21/04/2025 16:45
I set the vertex finding algorithm to create vertices from the front plane information with parameters (sigma = 0.5 ,n_iters = 25) to investigate why increasing the number of iterations was causing worse performance. In short, I still don't know, but here are some plots:
Basically the algorithm creates k vertices, it tests the range [1, num_endpoints] and does n_iters iterations and picks the one the one that minimizes the BIC score. I.e. it will still test every k even if the "minimizing" k is smaller than num_endpoints.
So most the time, the algorithm converges very quickly. But for some cases the convergence does not happen quickly. It may be that there is some very high order (insiginificant adjustment) that sometimes gets made at later iterations, but it's unclear from these plots.
21/04/2025 16:50
I wanted to get a deeper view of how the fit performance affected the pattern finding.
Note for all the below plots, we are just fitting the (x,z) or "front" plane. That may be made obvious in the second plot.
We define "failed" fits as just ones that did not have at least 2 unique z coordiantes to fit to. Unsurprisingly, the fewer failed fits there are, the better we do.
I also checked how the performance goes with number of muon hits. Also unsurprisingly, the more muon hits we see in the (x,z) plane, the better performance we have building patterns out of that (x,z) data.
21/04/2025 16:53
I wanted to try a new vertex seeding algorithm to see how it affects performance. I tried some bad ones (not shown) at first and saw a big drop in performance (~90% --> ~80% validation success rate for \boldsymbol{\pi^+ \rightarrow \mu^+ + \nu_\mu \rightarrow e^+ + \nu_e}).
My "new" algorithm basically does a round of clustering at first to determine vertex guesses. It just groups endpoints by distance, ensuring that two endpoints from the same tracklet are not grouped and assigns cluster centers based on that distance. The benefit to this is no random guess is made, the DBSCAN algorithm makes no such random guess in it's clustering. Any additional vertices requested (k are requested by the constrained_k_means algorithm) are just added randomly as before.
\boldsymbol{\pi^+ \rightarrow \mu^+ + \nu_\mu \rightarrow e^+ + \nu_e}:
\boldsymbol{\pi^+ \rightarrow e^+ + \nu_e}:
We see a performance increase in identifying "expected" events, but a drop-off for "unexpected" events like those with electrons.
22/04/2025 01:15
I pulled Jessie's tracklet finding algorithm and ran the simulation with this macro:
# print macro commands on screen /control/verbose 0 ################################################################### # geometry must be specified before /run/initialize # /geometry/source default_lyso.gdml # Configure magnetic fields # /magfield/scalingfactor 1 # End of geometry configuration # ################################################################### ################################################################### # Configuration of the physics to be used # /Physics/SelectList QGSP_BERT_EM4 #/Physics/SelectList CEX /process/had/verbose 0 /process/em/verbose 0 # Add optical physics. #/Physics/AddOptics #/process/optical/verbose 1 #/process/optical/processActivation Cerenkov false #/process/optical/processActivation Scintillation false #/process/optical/processActivation OpAbsorption true #/process/optical/processActivation OpRayleigh true #/process/optical/processActivation OpMieHG true #/process/optical/processActivation OpBoundary true #/process/optical/processActivation OpWLS true #/process/optical/processActivation OpWLS2 true # #/process/optical/scintillation/verbose 1 #/process/optical/scintillation/setByParticleType false #/process/optical/scintillation/setTrackInfo true #/process/optical/scintillation/setFiniteRiseTime false #/process/optical/scintillation/setStackPhotons true #/process/optical/scintillation/setTrackSecondariesFirst true # Decay mode selection #/decay/all /decay/pimunu #/decay/rad_muon #/decay/pienu #/decay/pienug #/decay/rad_michel #/decay/rad_michel_rad_muon #/decay/pibeta # Cap muon/pion lifetime to enhance probability of decay in flight # /decay/mulifetimecap 50 ps # /decay/pilifetimecap 100 ps # /decay/pistartatcreation false # Increase Pion Charge Exchange Cross Section by Factor #/Physics/scalePionCEX 1 # End of physics configuration # ################################################################### ################################################################### # Configuration of the output to be written # # path to output file. "#RUN#" will be replaced by the run ID /output/FileName ./pimunu_run#RUN# # Switching on/off branches in the output TTree /output/SaveInit true #/output/SaveTrack true /output/SaveDecay true /output/SaveAtar true /output/SaveTracker true /output/SaveCalo true /output/SaveSipm true /output/SaveUpstream true #/output/SaveGhost true #/output/SaveGhostCalo true #/output/SaveSplitoff true # End of output configuration # ################################################################### #================================================================== # Initialise the run /run/initialize # check physics processes and particles. # Beware, output is somewhat messy in multithreaded mode #/process/list #/particle/list #/geometry/list # Configure pion beam /gen/select beam # Select beam phasespace definition # Options: pre-built, flexible /gen/beam pre-built # Beam contaminations (0 - 1.00) /gen/beam/muons 0 /gen/beam/positrons 0 # Beam parameters # General beam parameters #/gen/beam/momentum 65 MeV #/gen/beam/momsigma 1.4 MeV #/gen/beam/xmean 0 mm #/gen/beam/ymean 0 mm #/gen/beam/zpos -1000 mm # Select gaussian beam and set mode # Options are (gaus, shape) for shape #/gen/beam/shape cyl # Cylindrical beam parameters #/gen/beam/cyl/rmax 10 mm # Gaussian beam parameters #/gen/beam/gaus/rmax 25 cm #/gen/beam/gaus/mode sinit_emittance #/gen/beam/gaus/zoff 0 mm #/gen/beam/gaus/xinitsigma 91 mm #/gen/beam/gaus/yinitsigma 54 mm #/gen/beam/gaus/xemittance 0.62 mm #/gen/beam/gaus/yemittance 0.23 mm #/gen/beam/gaus/xwaistsigma 6.8 mm #/gen/beam/gaus/ywaistsigma 4.3 mm #/gen/beam/gaus/xprimesigma 0.09 #/gen/beam/gaus/yprimesigma 0.05 #/gen/select signal #/gen/signal/momentum 70 MeV #/gen/signal/momsigma 1 MeV #/gen/signal/thetaMin 60 deg #/gen/signal/thetaMax 120 deg #/gen/signal/phiMin 0 deg #/gen/signal/phiMax 360 deg #/gen/signal/sigmaX 2 mm #/gen/signal/sigmaY 2 mm #/gen/signal/sigmaZ 2 mm # configure the generic particle source #/gps/particle pi+ #/gps/energy 0.0 MeV #/gps/pos/type Volume #/gps/pos/shape Para #/gps/pos/centre 0 0 0 mm #/gps/pos/halfx 10 mm #/gps/pos/halfy 10 mm #/gps/pos/halfz 3 mm # ================================================================= # visualize geometry and events for debugging #/vis/open HepRepFile #/vis/drawVolume #/vis/scene/add/trajectories # +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ # Start the run /run/beamOn 100000
# print macro commands on screen
/control/verbose 0
###################################################################
# geometry must be specified before /run/initialize #
/geometry/source default_lyso.gdml
# Configure magnetic fields
# /magfield/scalingfactor 1
# End of geometry configuration #
###################################################################
###################################################################
# Configuration of the physics to be used #
/Physics/SelectList QGSP_BERT_EM4
#/Physics/SelectList CEX
/process/had/verbose 0
/process/em/verbose 0
# Add optical physics.
#/Physics/AddOptics
#/process/optical/verbose 1
#/process/optical/processActivation Cerenkov false
#/process/optical/processActivation Scintillation false
#/process/optical/processActivation OpAbsorption true
#/process/optical/processActivation OpRayleigh true
#/process/optical/processActivation OpMieHG true
#/process/optical/processActivation OpBoundary true
#/process/optical/processActivation OpWLS true
#/process/optical/processActivation OpWLS2 true
#
#/process/optical/scintillation/verbose 1
#/process/optical/scintillation/setByParticleType false
#/process/optical/scintillation/setTrackInfo true
#/process/optical/scintillation/setFiniteRiseTime false
#/process/optical/scintillation/setStackPhotons true
#/process/optical/scintillation/setTrackSecondariesFirst true
# Decay mode selection
#/decay/all
/decay/pimunu
#/decay/rad_muon
#/decay/pienu
#/decay/pienug
#/decay/rad_michel
#/decay/rad_michel_rad_muon
#/decay/pibeta
# Cap muon/pion lifetime to enhance probability of decay in flight
# /decay/mulifetimecap 50 ps
# /decay/pilifetimecap 100 ps
# /decay/pistartatcreation false
# Increase Pion Charge Exchange Cross Section by Factor
#/Physics/scalePionCEX 1
# End of physics configuration #
###################################################################
###################################################################
# Configuration of the output to be written #
# path to output file. "#RUN#" will be replaced by the run ID
/output/FileName ./pimunu_run#RUN#
# Switching on/off branches in the output TTree
/output/SaveInit true
#/output/SaveTrack true
/output/SaveDecay true
/output/SaveAtar true
/output/SaveTracker true
/output/SaveCalo true
/output/SaveSipm true
/output/SaveUpstream true
#/output/SaveGhost true
#/output/SaveGhostCalo true
#/output/SaveSplitoff true
# End of output configuration #
###################################################################
#==================================================================
# Initialise the run
/run/initialize
# check physics processes and particles.
# Beware, output is somewhat messy in multithreaded mode
#/process/list
#/particle/list
#/geometry/list
# Configure pion beam
/gen/select beam
# Select beam phasespace definition
# Options: pre-built, flexible
/gen/beam pre-built
# Beam contaminations (0 - 1.00)
/gen/beam/muons 0
/gen/beam/positrons 0
# Beam parameters
# General beam parameters
#/gen/beam/momentum 65 MeV
#/gen/beam/momsigma 1.4 MeV
#/gen/beam/xmean 0 mm
#/gen/beam/ymean 0 mm
#/gen/beam/zpos -1000 mm
# Select gaussian beam and set mode
# Options are (gaus, shape) for shape
#/gen/beam/shape cyl
# Cylindrical beam parameters
#/gen/beam/cyl/rmax 10 mm
# Gaussian beam parameters
#/gen/beam/gaus/rmax 25 cm
#/gen/beam/gaus/mode sinit_emittance
#/gen/beam/gaus/zoff 0 mm
#/gen/beam/gaus/xinitsigma 91 mm
#/gen/beam/gaus/yinitsigma 54 mm
#/gen/beam/gaus/xemittance 0.62 mm
#/gen/beam/gaus/yemittance 0.23 mm
#/gen/beam/gaus/xwaistsigma 6.8 mm
#/gen/beam/gaus/ywaistsigma 4.3 mm
#/gen/beam/gaus/xprimesigma 0.09
#/gen/beam/gaus/yprimesigma 0.05
#/gen/select signal
#/gen/signal/momentum 70 MeV
#/gen/signal/momsigma 1 MeV
#/gen/signal/thetaMin 60 deg
#/gen/signal/thetaMax 120 deg
#/gen/signal/phiMin 0 deg
#/gen/signal/phiMax 360 deg
#/gen/signal/sigmaX 2 mm
#/gen/signal/sigmaY 2 mm
#/gen/signal/sigmaZ 2 mm
# configure the generic particle source
#/gps/particle pi+
#/gps/energy 0.0 MeV
#/gps/pos/type Volume
#/gps/pos/shape Para
#/gps/pos/centre 0 0 0 mm
#/gps/pos/halfx 10 mm
#/gps/pos/halfy 10 mm
#/gps/pos/halfz 3 mm
# =================================================================
# visualize geometry and events for debugging
#/vis/open HepRepFile
#/vis/drawVolume
#/vis/scene/add/trajectories
# +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# Start the run
/run/beamOn 100000
As shown below, it appears there are fewer events with lots of "extra" particles. Maybe this is just chance, maybe I've somehow repressed them with these new settings.
22/04/2025 01:05
Some results using Jessie's endpoints; However, this is just using the x-z information and random tracklet seeding (also some endpoints can be NaN for \boldsymbol{\pi^+ \rightarrow \mu^+ + \nu_\mu \rightarrow e^+ + \nu_e} decays :
Instead using 3D information: (I still have the "worst of both worlds" theory here, where you're now twice as likely to have an endpoint have a messed up coordinate, causing issues)
Now if we use front information and a non-random seeding:
And 3D information with non-random seeding
22/04/2025 01:08
Here's the same plots for \boldsymbol{\pi^+ \rightarrow e^+ + \nu_e}.
Again, the first set is just (x,z) information and random tracklet seeding:
If we instead use a 3D information (also some endpoints can be NaN):
Now if we use 3D information a non-random seeding (remapped NaN --> (0,0,0) as well):
22/04/2025 14:45
For reference, here are the reuslts using the same dataset, but the old endpoint finding algorithm (only using front planes, \boldsymbol{\pi^+ \rightarrow \mu^+ + \nu_\mu \rightarrow e^+ + \nu_e}
So we see \frac{98.65}{92.57} \approx 1.07 \implies 7% increase in performance.
What's not shown is that these new endpoints also solve the "missing tracklet problems" I described earlier (you can sort of see this affect by looking that the pattern correct is now never true while validation is false). This indicates the problem was with my endpoint assigning (probably because I had to assign some endpoints to "None" and the vertex finder didn't like that).
23/04/2025 15:07
All of the above data was taken with respect to the truth pattern finder running on the reconstructed tracklets This makes the pattern finding performance look as if it is better than it truly is. It also makes the the reconstructed pattern particle compoisition look the same as the true particle composition, because the pattern finder only has access to the reco tracklets knowledge of particle compositon.
Instead if we comapre to truth pattern finding running on truth tracklets, we get different results:
\boldsymbol{\pi^+ \rightarrow \mu^+ + \nu_\mu \rightarrow e^+ + \nu_e}
(x,z) information, random seeding, Jessie's Tracklets:
(x,z) information, random seeding, fitting truth Tracklets:
In a way, this differentiates how the pattern finding algorithm and tracklet finding algoirthm affect performance. In any event, you still see performance gains using Jessie's Tracklet finding, which is good.
23/04/2025 15:36
Here's the same plots for \boldsymbol{\pi^+ \rightarrow e^+ + \nu_e}.
(x,z) information, random seeding, Jessie's Tracklets:
(x,z) information, random seeding, fitting truth Tracklets:
23/04/2025 17:18
It seems some endpoints just don't have enough information to form at all (from inspection, I don't see cases where where have an (x,z) coordinate but not (y,z) coordiante)
⚠️ Invalid endpoint detected: [ nan nan -inf] [array([ 4.98903478, -6.37736092, 0.24183216])] ⚠️ Invalid endpoint at ([ nan nan -inf]) for vertex 0. ⚠️ Non-finite values detected in the input endpoints and cluster centers: k: 1 End Points: [[ 5.02008748 -6.20570401 -0.15888712] [ 4.69999994 -6.91760125 3.42886032] [ nan nan -inf] [ 4.69999993 -6.86386318 3.41874448]] Cluster Centers: [[ 4.80669578 -6.66238948 2.22957256] [ 4.80669578 -6.66238948 2.22957256] [ 4.80669578 -6.66238948 2.22957256] [ 4.80669578 -6.66238948 2.22957256]] ⚠️ Non-finite values detected in a vector pair! Cluster Centers: [ 4.80669578 -6.66238948 2.22957256] End Points: [ nan nan -inf] Valid mask: [False False False] Filtered vec1: [] Filtered vec2: []
⚠️ Invalid endpoint detected: [ nan nan -inf]
[array([ 4.98903478, -6.37736092, 0.24183216])]
⚠️ Invalid endpoint at ([ nan nan -inf]) for vertex 0.
⚠️ Non-finite values detected in the input endpoints and cluster centers:
k: 1
End Points: [[ 5.02008748 -6.20570401 -0.15888712]
[ 4.69999994 -6.91760125 3.42886032]
[ nan nan -inf]
[ 4.69999993 -6.86386318 3.41874448]]
Cluster Centers: [[ 4.80669578 -6.66238948 2.22957256]
[ 4.80669578 -6.66238948 2.22957256]
[ 4.80669578 -6.66238948 2.22957256]
[ 4.80669578 -6.66238948 2.22957256]]
⚠️ Non-finite values detected in a vector pair!
Cluster Centers: [ 4.80669578 -6.66238948 2.22957256]
End Points: [ nan nan -inf]
Valid mask: [False False False]
Filtered vec1: []
Filtered vec2: []I added some code to "ignore" these endpoints when creating clusters, but it doesn't seem to imrpove anything. I.e. the code as it was before was already ignoring these endpoints because they produced an infinite BIC.
Performance for \boldsymbol{\pi^+ \rightarrow e^+ + \nu_e} with endpoint filtering. Using Jessie's Tracklets, Both planes of information, random seeding:
Performance for \boldsymbol{\pi^+ \rightarrow \mu^+ + \nu_\mu \rightarrow e^+ + \nu_e} with endpoint filtering. Using Jessie's Tracklets, Both planes of information, random seeding:
Currently erroring out because some points get deleted and then later attempted to access